﻿/******************************************************************************
 * 
 * Announce: CSharpKit, Basic algorithms, components and definitions.
 *           Copyright (C) ShenYongchen.
 *           All rights reserved.
 *   Author: 申永辰.郑州 (shenyczz@163.com)
 *  WebSite: http://github.com/shenyczz/CSharpKit
 *
 * THIS CODE IS LICENSED UNDER THE MIT LICENSE (MIT).
 * THIS CODE IS PROVIDED *AS IS* WITHOUT WARRANTY OF 
 * ANY KIND, EITHER EXPRESS OR IMPLIED, INCLUDING ANY
 * IMPLIED WARRANTIES OF FITNESS FOR A PARTICULAR
 * PURPOSE, MERCHANTABILITY, OR NON-INFRINGEMENT.
 * 
******************************************************************************/

using System;
using System.Collections.Generic;
using System.Linq;
using System.Runtime.Serialization;
using CSharpKit.Threading.Tasks;

namespace CSharpKit.Numerics.LinearAlgebra
{
    public interface IMatrixStorage : IStorage { }

    /// <summary>
    /// 矩阵存储
    /// </summary>
    /// <typeparam name="T"></typeparam>
    [Serializable]
    public partial class MatrixStorage<T> : Numeric<T>, IMatrixStorage, IEquatable<MatrixStorage<T>>
        where T : struct, IFormattable, IEquatable<T>
    {
        #region Constructors

        protected MatrixStorage()
            : this(0, 0) { }

        public MatrixStorage(int row, int column)
        {
            if (row <= 0)
            {
                throw new ArgumentOutOfRangeException(nameof(row), Resources.Numerics.MatrixRowsMustBePositive);
            }

            if (column <= 0)
            {
                throw new ArgumentOutOfRangeException(nameof(column), Resources.Numerics.MatrixColumnsMustBePositive);
            }


            _Rows = row;
            _Cols = column;

            _IsDense = true;
            _IsMutableFully = true;
            _Order = Matrix<T>.MajorOrder;
        }

        #endregion


        #region Cols

        readonly int _Cols;
        [DataMember(Order = 2)]
        public int ColumnCount => _Cols;

        #endregion

        #region Rows

        readonly int _Rows;
        [DataMember(Order = 1)]
        public int RowCount => _Rows;

        #endregion

        #region Order

        protected MajorOrders _Order;
        public MajorOrders Order => _Order;

        #endregion

        #region IsDense

        protected bool _IsDense;
        /// <summary>
        /// 如果矩阵存储格式密集，则为真。
        /// </summary>
        public bool IsDense => _IsDense;

        #endregion

        #region IsMutableFully

        protected bool _IsMutableFully;
        /// <summary>
        /// 全部易变<br/>
        /// 矩阵的所有字段都可以设置为任何值为 True <br/>
        /// 有些字段是固定的，就像在对角线矩阵上一样则为 False
        /// </summary>
        public bool IsMutableFully => _IsMutableFully;

        #endregion

        #region IsMutableAt

        public virtual bool IsMutableAt(int row, int column) => true;

        #endregion


        #region Indexer - 索引器

        public T this[int row, int column]
        {
            get
            {
                ValidateRange(row, column);
                return At(row, column);
            }
            set
            {
                ValidateRange(row, column);
                At(row, column, value);
            }
        }

        internal virtual T At(int row, int column) => throw new Exception();
        internal virtual void At(int row, int column, T tvalue) => throw new Exception();

        /// <summary>
        /// 元素位置 <br/>
        /// 行主序: pos = (row * Cols) + col <br/>
        /// 列主序: pos = (col * Rows) + row <br/>
        /// </summary>
        /// <param name="row"></param>
        /// <param name="col"></param>
        /// <returns></returns>
        internal virtual int ElementPosition(int row, int col) => throw new Exception();
        /// <summary>
        /// Evaluate the row and column at a specific data index.
        /// 
        /// 行主序
        /// row = index / ColumnCount;
        /// col = index % ColumnCount;
        /// row = Math.DivRem(index, ColumnCount, out col);
        /// 
        ///  列主序
        /// row = index % RowCount;
        /// col = index / RowCount;
        /// column = Math.DivRem(index, RowCount, out row);
        /// </summary>
        protected virtual void RowColumnAtIndex(int index, out int row, out int col) => throw new Exception();

        #endregion


        #region IEquatable<MatrixStorage<T>>

        public override bool Equals(object obj)
        {
            return Equals(obj as MatrixStorage<T>);
        }

        public bool Equals(MatrixStorage<T> other)
        {
            if (other == null)
            {
                return false;
            }

            if (RowCount != other.RowCount)
            {
                return false;
            }

            if (ColumnCount != other.ColumnCount)
            {
                return false;
            }

            // Accept if the argument is the same object as this.
            if (ReferenceEquals(this, other))
            {
                return true;
            }

            // Perform element wise comparison.
            return Find2Unchecked(other, (a, b) => !a.Equals(b), Zeros.AllowSkip) == null;
        }

        public override int GetHashCode()
        {
            return base.GetHashCode()
                ^ RowCount.GetHashCode()
                ^ ColumnCount.GetHashCode()
                ^ TheGetHashCode()
                ;
        }
        int TheGetHashCode()
        {
            var hashNum = System.Math.Min(RowCount * ColumnCount, 25);
            int hash = 17;
            unchecked
            {
                for (var i = 0; i < hashNum; i++)
                {
                    int col = i % ColumnCount;
                    int row = i / ColumnCount;
                    //int row = Math.DivRem(i, ColumnCount, out int col);
                    hash = hash * 31 + At(row, col).GetHashCode();
                }
            }
            return hash;
        }

        #endregion


        #region Clearing

        public virtual void Clear()
        {
            for (var i = 0; i < RowCount; i++)
            {
                for (var j = 0; j < ColumnCount; j++)
                {
                    At(i, j, Zero);
                }
            }
        }
        public void Clear(int rowIndex, int rowCount, int columnIndex, int columnCount)
        {
            if (rowCount < 1 || columnCount < 1)
            {
                return;
            }

            if (rowIndex + rowCount > RowCount || rowIndex < 0)
            {
                throw new ArgumentOutOfRangeException(nameof(rowIndex));
            }

            if (columnIndex + columnCount > ColumnCount || columnIndex < 0)
            {
                throw new ArgumentOutOfRangeException(nameof(columnIndex));
            }

            ClearUnchecked(rowIndex, rowCount, columnIndex, columnCount);
        }
        public void ClearRows(int[] rowIndices)
        {
            if (rowIndices.Length == 0)
            {
                return;
            }

            for (int k = 0; k < rowIndices.Length; k++)
            {
                if (rowIndices[k] < 0 || rowIndices[k] >= RowCount)
                {
                    throw new ArgumentOutOfRangeException(nameof(rowIndices));
                }
            }

            ClearRowsUnchecked(rowIndices);
        }
        public void ClearColumns(int[] columnIndices)
        {
            if (columnIndices.Length == 0)
            {
                return;
            }

            for (int k = 0; k < columnIndices.Length; k++)
            {
                if ((uint)columnIndices[k] >= (uint)ColumnCount)
                {
                    throw new ArgumentOutOfRangeException(nameof(columnIndices));
                }
            }

            ClearColumnsUnchecked(columnIndices);
        }

        /// <summary>
        /// 清除（无错误检查）
        /// 调用前需要进行数组边界检查
        /// </summary>
        /// <param name="rowIndex"></param>
        /// <param name="rowCount"></param>
        /// <param name="columnIndex"></param>
        /// <param name="columnCount"></param>
        internal virtual void ClearUnchecked(int rowIndex, int rowCount, int columnIndex, int columnCount)
        {
            for (var i = rowIndex; i < rowIndex + rowCount; i++)
            {
                for (var j = columnIndex; j < columnIndex + columnCount; j++)
                {
                    At(i, j, Zero);
                }
            }
        }
        /// <summary>
        /// 清除行元素（无错误检查）
        /// </summary>
        /// <param name="rowIndices">行索引数组</param>
        internal virtual void ClearRowsUnchecked(int[] rowIndices)
        {
            for (var k = 0; k < rowIndices.Length; k++)
            {
                int row = rowIndices[k];
                for (var j = 0; j < ColumnCount; j++)
                {
                    At(row, j, Zero);
                }
            }
        }
        /// <summary>
        /// 清除列元素（无错误检查）
        /// </summary>
        /// <param name="columnIndices">列索引数组</param>
        internal virtual void ClearColumnsUnchecked(int[] columnIndices)
        {
            for (var k = 0; k < columnIndices.Length; k++)
            {
                int column = columnIndices[k];
                for (var i = 0; i < RowCount; i++)
                {
                    At(i, column, Zero);
                }
            }
        }

        #endregion


        #region Copying

        // MATRIX COPY

        public void CopyTo(MatrixStorage<T> target, ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (ReferenceEquals(this, target))
            {
                return;
            }

            if (RowCount != target.RowCount || ColumnCount != target.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, target.RowCount + "x" + target.ColumnCount);
                throw new ArgumentException(message, nameof(target));
            }

            CopyToUnchecked(target, existingData);
        }

        internal virtual void CopyToUnchecked(MatrixStorage<T> target, ExistingData existingData)
        {
            for (int j = 0; j < ColumnCount; j++)
            {
                for (int i = 0; i < RowCount; i++)
                {
                    target.At(i, j, At(i, j));
                }
            }
        }

        public void CopySubMatrixTo(MatrixStorage<T> target,
            int sourceRowIndex, int targetRowIndex, int rowCount,
            int sourceColumnIndex, int targetColumnIndex, int columnCount,
            ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (rowCount == 0 || columnCount == 0)
            {
                return;
            }

            if (ReferenceEquals(this, target))
            {
                throw new NotSupportedException();
            }

            ValidateSubMatrixRange(target,
                sourceRowIndex, targetRowIndex, rowCount,
                sourceColumnIndex, targetColumnIndex, columnCount);

            CopySubMatrixToUnchecked(target, sourceRowIndex, targetRowIndex, rowCount,
                sourceColumnIndex, targetColumnIndex, columnCount, existingData);
        }

        internal virtual void CopySubMatrixToUnchecked(MatrixStorage<T> target,
            int sourceRowIndex, int targetRowIndex, int rowCount,
            int sourceColumnIndex, int targetColumnIndex, int columnCount,
            ExistingData existingData)
        {
            for (int j = sourceColumnIndex, jj = targetColumnIndex; j < sourceColumnIndex + columnCount; j++, jj++)
            {
                for (int i = sourceRowIndex, ii = targetRowIndex; i < sourceRowIndex + rowCount; i++, ii++)
                {
                    target.At(ii, jj, At(i, j));
                }
            }
        }

        // ROW COPY

        public void CopyRowTo(VectorStorage<T> target, int rowIndex, ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            ValidateRowRange(target, rowIndex);
            CopySubRowToUnchecked(target, rowIndex, 0, 0, ColumnCount, existingData);
        }
        internal virtual void CopySubRowToUnchecked(VectorStorage<T> target, int rowIndex,
            int sourceColumnIndex, int targetColumnIndex, int columnCount, ExistingData existingData)
        {
            for (int j = sourceColumnIndex, jj = targetColumnIndex; j < sourceColumnIndex + columnCount; j++, jj++)
            {
                target.At(jj, At(rowIndex, j));
            }
        }

        public void CopySubRowTo(VectorStorage<T> target, int rowIndex,
            int sourceColumnIndex, int targetColumnIndex, int columnCount,
            ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (columnCount == 0)
            {
                return;
            }

            ValidateSubRowRange(target, rowIndex, sourceColumnIndex, targetColumnIndex, columnCount);
            CopySubRowToUnchecked(target, rowIndex, sourceColumnIndex, targetColumnIndex, columnCount, existingData);
        }

        // COLUMN COPY

        public void CopyColumnTo(VectorStorage<T> target, int columnIndex, ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            ValidateColumnRange(target, columnIndex);
            CopySubColumnToUnchecked(target, columnIndex, 0, 0, RowCount, existingData);
        }

        public void CopySubColumnTo(VectorStorage<T> target, int columnIndex,
            int sourceRowIndex, int targetRowIndex, int rowCount,
            ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (rowCount == 0)
            {
                return;
            }

            ValidateSubColumnRange(target, columnIndex, sourceRowIndex, targetRowIndex, rowCount);
            CopySubColumnToUnchecked(target, columnIndex, sourceRowIndex, targetRowIndex, rowCount, existingData);
        }

        internal virtual void CopySubColumnToUnchecked(VectorStorage<T> target, int columnIndex,
            int sourceRowIndex, int targetRowIndex, int rowCount, ExistingData existingData)
        {
            for (int i = sourceRowIndex, ii = targetRowIndex; i < sourceRowIndex + rowCount; i++, ii++)
            {
                target.At(ii, At(i, columnIndex));
            }
        }



        #endregion


        #region Enumerating

        public virtual IEnumerable<T> Enumerate()
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    yield return At(i, j);
                }
            }
        }
        public virtual IEnumerable<T> EnumerateNonZero()
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    var x = At(i, j);
                    if (!Zero.Equals(x))
                    {
                        yield return x;
                    }
                }
            }
        }
        public virtual IEnumerable<Tuple<int, int, T>> EnumerateIndexed()
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    yield return new Tuple<int, int, T>(i, j, At(i, j));
                }
            }
        }
        public virtual IEnumerable<Tuple<int, int, T>> EnumerateNonZeroIndexed()
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    var x = At(i, j);
                    if (!Zero.Equals(x))
                    {
                        yield return new Tuple<int, int, T>(i, j, x);
                    }
                }
            }
        }

        #endregion


        #region Extracting

        public virtual T[,] ToArray()
        {
            var result = new T[RowCount, ColumnCount];

            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    result[i, j] = At(i, j);
                }
            }

            return result;
        }
        public virtual T[][] ToRowArrays()
        {
            var ret = new T[RowCount][];

            for (int i = 0; i < RowCount; i++)
            {
                var row = new T[ColumnCount];
                for (int j = 0; j < ColumnCount; j++)
                {
                    row[j] = At(i, j);
                }
                ret[i] = row;
            }

            return ret;
        }
        public virtual T[][] ToColumnArrays()
        {
            var ret = new T[ColumnCount][];
            for (int j = 0; j < ColumnCount; j++)
            {
                var column = new T[RowCount];
                for (int i = 0; i < RowCount; i++)
                {
                    column[i] = At(i, j);
                }
                ret[j] = column;
            }
            return ret;
        }
        public virtual T[] ToRowMajorArray()
        {
            var ret = new T[RowCount * ColumnCount];
            for (int i = 0; i < RowCount; i++)
            {
                var offset = i * ColumnCount;
                for (int j = 0; j < ColumnCount; j++)
                {
                    ret[offset + j] = At(i, j);
                }
            }
            return ret;
        }
        public virtual T[] ToColumnMajorArray()
        {
            var ret = new T[RowCount * ColumnCount];
            for (int j = 0; j < ColumnCount; j++)
            {
                var offset = j * RowCount;
                for (int i = 0; i < RowCount; i++)
                {
                    ret[offset + i] = At(i, j);
                }
            }
            return ret;
        }


        #endregion


        #region Finding

        public virtual Tuple<int, int, T> Find(Func<T, bool> predicate, Zeros zeros)
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    var item = At(i, j);
                    if (predicate(item))
                    {
                        return new Tuple<int, int, T>(i, j, item);
                    }
                }
            }
            return null;
        }
        public Tuple<int, int, T, TOther> Find2<TOther>(MatrixStorage<TOther> other, Func<T, TOther, bool> predicate, Zeros zeros)
            where TOther : struct, IEquatable<TOther>, IFormattable
        {
            if (other == null)
            {
                throw new ArgumentNullException(nameof(other));
            }

            if (RowCount != other.RowCount || ColumnCount != other.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, other.RowCount + "x" + other.ColumnCount);
                throw new ArgumentException(message, nameof(other));
            }

            return Find2Unchecked(other, predicate, zeros);
        }
        internal virtual Tuple<int, int, T, TOther> Find2Unchecked<TOther>(MatrixStorage<TOther> other, Func<T, TOther, bool> predicate, Zeros zeros)
            where TOther : struct, IEquatable<TOther>, IFormattable
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    var item = At(i, j);
                    var otherItem = other.At(i, j);
                    if (predicate(item, otherItem))
                    {
                        return new Tuple<int, int, T, TOther>(i, j, item, otherItem);
                    }
                }
            }

            return null;
        }

        #endregion


        #region COMBINATORS - Maping

        public virtual void MapInplace(Func<T, T> fun, Zeros zeros)
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    At(i, j, fun(At(i, j)));
                }
            }
        }

        public virtual void MapIndexedInplace(Func<int, int, T, T> f, Zeros zeros)
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    At(i, j, f(i, j, At(i, j)));
                }
            }
        }

        public virtual void MapTo<TU>(MatrixStorage<TU> target, Func<T, TU> f, Zeros zeros, ExistingData existingData)
            where TU : struct, IFormattable, IEquatable<TU>
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (RowCount != target.RowCount || ColumnCount != target.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, target.RowCount + "x" + target.ColumnCount);
                throw new ArgumentException(message, nameof(target));
            }

            MapToUnchecked(target, f, zeros, existingData);
        }

        public virtual void Map2To(MatrixStorage<T> target, MatrixStorage<T> other, Func<T, T, T> fun, Zeros zeros, ExistingData existingData)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (other == null)
            {
                throw new ArgumentNullException(nameof(other));
            }

            if (RowCount != target.RowCount || ColumnCount != target.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, target.RowCount + "x" + target.ColumnCount);
                throw new ArgumentException(message, nameof(target));
            }

            if (RowCount != other.RowCount || ColumnCount != other.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, other.RowCount + "x" + other.ColumnCount);
                throw new ArgumentException(message, nameof(other));
            }

            Map2ToUnchecked(target, other, fun, zeros, existingData);
        }

        public void MapIndexedTo<TU>(MatrixStorage<TU> target, Func<int, int, T, TU> f, Zeros zeros, ExistingData existingData)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (RowCount != target.RowCount || ColumnCount != target.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, target.RowCount + "x" + target.ColumnCount);
                throw new ArgumentException(message, nameof(target));
            }

            MapIndexedToUnchecked(target, f, zeros, existingData);
        }

        public void MapSubMatrixIndexedTo<TU>(MatrixStorage<TU> target, Func<int, int, T, TU> f,
            int sourceRowIndex, int targetRowIndex, int rowCount,
            int sourceColumnIndex, int targetColumnIndex, int columnCount,
            Zeros zeros, ExistingData existingData)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (rowCount == 0 || columnCount == 0)
            {
                return;
            }

            if (ReferenceEquals(this, target))
            {
                throw new NotSupportedException();
            }

            ValidateSubMatrixRange(target,
                sourceRowIndex, targetRowIndex, rowCount,
                sourceColumnIndex, targetColumnIndex, columnCount);

            MapSubMatrixIndexedToUnchecked(target, f, sourceRowIndex, targetRowIndex, rowCount, sourceColumnIndex, targetColumnIndex, columnCount, zeros, existingData);
        }

        internal virtual void MapToUnchecked<TU>(MatrixStorage<TU> target, Func<T, TU> f, Zeros zeros, ExistingData existingData)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    target.At(i, j, f(At(i, j)));
                }
            }
        }
        internal virtual void Map2ToUnchecked(MatrixStorage<T> target, MatrixStorage<T> other, Func<T, T, T> fun, Zeros zeros, ExistingData existingData)
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    target.At(i, j, fun(At(i, j), other.At(i, j)));
                }
            }
        }
        internal virtual void MapIndexedToUnchecked<TU>(MatrixStorage<TU> target, Func<int, int, T, TU> f, Zeros zeros, ExistingData existingData)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            for (int j = 0; j < ColumnCount; j++)
            {
                for (int i = 0; i < RowCount; i++)
                {
                    target.At(i, j, f(i, j, At(i, j)));
                }
            }
        }
        internal virtual void MapSubMatrixIndexedToUnchecked<TU>(MatrixStorage<TU> target, Func<int, int, T, TU> f,
            int sourceRowIndex, int targetRowIndex, int rowCount,
            int sourceColumnIndex, int targetColumnIndex, int columnCount,
            Zeros zeros, ExistingData existingData)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            for (int j = sourceColumnIndex, jj = targetColumnIndex; j < sourceColumnIndex + columnCount; j++, jj++)
            {
                for (int i = sourceRowIndex, ii = targetRowIndex; i < sourceRowIndex + rowCount; i++, ii++)
                {
                    target.At(ii, jj, f(ii, jj, At(i, j)));
                }
            }
        }


        #endregion


        #region COMBINATORS - Folding

        public void FoldByRow<TU>(TU[] target, Func<TU, T, TU> f, Func<TU, int, TU> finalize, TU[] state, Zeros zeros)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }
            if (target.Length != RowCount)
            {
                throw new ArgumentException(Resources.Numerics.ArgumentVectorsSameLength, nameof(target));
            }

            if (state == null)
            {
                throw new ArgumentNullException(nameof(state));
            }
            if (state.Length != RowCount)
            {
                throw new ArgumentException(Resources.Numerics.ArgumentVectorsSameLength, nameof(state));
            }

            FoldByRowUnchecked(target, f, finalize, state, zeros);
        }
        public void FoldByColumn<TU>(TU[] target, Func<TU, T, TU> f, Func<TU, int, TU> finalize, TU[] state, Zeros zeros)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }
            if (target.Length != ColumnCount)
            {
                throw new ArgumentException(Resources.Numerics.ArgumentVectorsSameLength, nameof(target));
            }

            if (state == null)
            {
                throw new ArgumentNullException(nameof(state));
            }
            if (state.Length != ColumnCount)
            {
                throw new ArgumentException(Resources.Numerics.ArgumentVectorsSameLength, nameof(state));
            }

            FoldByColumnUnchecked(target, f, finalize, state, zeros);
        }
        public TState Fold2<TOther, TState>(MatrixStorage<TOther> other, Func<TState, T, TOther, TState> f, TState state, Zeros zeros)
            where TOther : struct, IEquatable<TOther>, IFormattable
        {
            if (other == null)
            {
                throw new ArgumentNullException(nameof(other));
            }

            if (RowCount != other.RowCount || ColumnCount != other.ColumnCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, other.RowCount + "x" + other.ColumnCount);
                throw new ArgumentException(message, nameof(other));
            }

            return Fold2Unchecked(other, f, state, zeros);
        }

        internal virtual void FoldByRowUnchecked<TU>(TU[] target, Func<TU, T, TU> f, Func<TU, int, TU> finalize, TU[] state, Zeros zeros)
        {
            for (int i = 0; i < RowCount; i++)
            {
                TU s = state[i];
                for (int j = 0; j < ColumnCount; j++)
                {
                    s = f(s, At(i, j));
                }
                target[i] = finalize(s, ColumnCount);
            }
        }
        internal virtual void FoldByColumnUnchecked<TU>(TU[] target, Func<TU, T, TU> f, Func<TU, int, TU> finalize, TU[] state, Zeros zeros)
        {
            for (int j = 0; j < ColumnCount; j++)
            {
                TU s = state[j];
                for (int i = 0; i < RowCount; i++)
                {
                    s = f(s, At(i, j));
                }
                target[j] = finalize(s, RowCount);
            }
        }
        internal virtual TState Fold2Unchecked<TOther, TState>(MatrixStorage<TOther> other, Func<TState, T, TOther, TState> f, TState state, Zeros zeros)
            where TOther : struct, IEquatable<TOther>, IFormattable
        {
            for (int i = 0; i < RowCount; i++)
            {
                for (int j = 0; j < ColumnCount; j++)
                {
                    state = f(state, At(i, j), other.At(i, j));
                }
            }

            return state;
        }

        #endregion



        #region Transpose - 转置

        public void TransposeTo(MatrixStorage<T> target, ExistingData existingData = ExistingData.Clear)
        {
            if (target == null)
            {
                throw new ArgumentNullException(nameof(target));
            }

            if (RowCount != target.ColumnCount || ColumnCount != target.RowCount)
            {
                var message = string.Format(Resources.Numerics.ArgumentMatrixDimensions2, RowCount + "x" + ColumnCount, target.RowCount + "x" + target.ColumnCount);
                throw new ArgumentException(message, nameof(target));
            }

            if (ReferenceEquals(this, target))
            {
                TransposeSquareInplaceUnchecked();
                return;
            }

            TransposeToUnchecked(target, existingData);
        }

        internal virtual void TransposeToUnchecked(MatrixStorage<T> target, ExistingData existingData)
        {
            for (int j = 0; j < ColumnCount; j++)
            {
                for (int i = 0; i < RowCount; i++)
                {
                    target.At(j, i, At(i, j));
                }
            }
        }
        internal virtual void TransposeSquareInplaceUnchecked()
        {
            for (int j = 0; j < ColumnCount; j++)
            {
                for (int i = 0; i < j; i++)
                {
                    T swap = At(i, j);
                    At(i, j, At(j, i));
                    At(j, i, swap);
                }
            }
        }

        #endregion


        public static MatrixStorage<T> OfValue(int rows, int columns, T value)
        {
            var storage = default(MatrixStorage<T>);

            var order = Matrix<T>.MajorOrder;

            switch (order)
            {
                case MajorOrders.Column:
                case MajorOrders.Row:
                    storage = new MatrixStorageDense<T>(rows, columns);
                    break;

                case MajorOrders.Diagonal:
                    storage = new MatrixStorageDiagonal<T>(rows, columns);
                    break;

                case MajorOrders.None:
                    break;

                default:
                    break;
            }

            // Enumerable.All()


            if (storage != null)
            {
                var data = (storage as MatrixStorageDense<T>).Data;
                CommonParallel.For(0, data.Length, 4096, (a, b) =>
                {
                    for (int i = a; i < b; i++)
                    {
                        data[i] = value;
                    }
                });
            }

            return storage;
        }







        #region Validation

        private void ValidateRange(int row, int column)
        {
            if ((uint)row >= (uint)RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(row));
            }

            if ((uint)column >= (uint)ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(column));
            }
        }

        private void ValidateRowRange<TU>(VectorStorage<TU> target, int rowIndex)
            where TU : struct, IFormattable, IEquatable<TU>
        {
            if ((uint)rowIndex >= (uint)RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(rowIndex));
            }

            if (ColumnCount != target.Length)
            {
                throw new ArgumentException(Resources.Numerics.ArgumentMatrixSameRowDimension, nameof(target));
            }
        }

        private void ValidateColumnRange<TU>(VectorStorage<TU> target, int columnIndex)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            if ((uint)columnIndex >= (uint)ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(columnIndex));
            }

            if (RowCount != target.Length)
            {
                throw new ArgumentException(Resources.Numerics.ArgumentMatrixSameColumnDimension, nameof(target));
            }
        }

        private void ValidateSubMatrixRange<TU>(MatrixStorage<TU> target,
            int sourceRowIndex, int targetRowIndex, int rowCount,
            int sourceColumnIndex, int targetColumnIndex, int columnCount)
            where TU : struct, IFormattable, IEquatable<TU>
        {
            if (rowCount < 1)
            {
                throw new ArgumentOutOfRangeException(nameof(rowCount), Resources.Numerics.ArgumentMustBePositive);
            }

            if (columnCount < 1)
            {
                throw new ArgumentOutOfRangeException(nameof(columnCount), Resources.Numerics.ArgumentMustBePositive);
            }

            // Verify Source

            if ((uint)sourceRowIndex >= (uint)RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(sourceRowIndex));
            }

            if ((uint)sourceColumnIndex >= (uint)ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(sourceColumnIndex));
            }

            var sourceRowMax = sourceRowIndex + rowCount;
            var sourceColumnMax = sourceColumnIndex + columnCount;

            if (sourceRowMax > RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(rowCount));
            }

            if (sourceColumnMax > ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(columnCount));
            }

            // Verify Target

            if ((uint)targetRowIndex >= (uint)target.RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(targetRowIndex));
            }

            if ((uint)targetColumnIndex >= (uint)target.ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(targetColumnIndex));
            }

            var targetRowMax = targetRowIndex + rowCount;
            var targetColumnMax = targetColumnIndex + columnCount;

            if (targetRowMax > target.RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(rowCount));
            }

            if (targetColumnMax > target.ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(columnCount));
            }
        }

        private void ValidateSubRowRange<TU>(VectorStorage<TU> target, int rowIndex,
            int sourceColumnIndex, int targetColumnIndex, int columnCount)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            if (columnCount < 1)
            {
                throw new ArgumentOutOfRangeException(nameof(columnCount), Resources.Numerics.ArgumentMustBePositive);
            }

            // Verify Source

            if ((uint)rowIndex >= (uint)RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(rowIndex));
            }

            if ((uint)sourceColumnIndex >= (uint)ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(sourceColumnIndex));
            }

            if (sourceColumnIndex + columnCount > ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(columnCount));
            }

            // Verify Target

            if ((uint)targetColumnIndex >= (uint)target.Length)
            {
                throw new ArgumentOutOfRangeException(nameof(targetColumnIndex));
            }

            if (targetColumnIndex + columnCount > target.Length)
            {
                throw new ArgumentOutOfRangeException(nameof(columnCount));
            }
        }

        private void ValidateSubColumnRange<TU>(VectorStorage<TU> target, int columnIndex,
            int sourceRowIndex, int targetRowIndex, int rowCount)
            where TU : struct, IEquatable<TU>, IFormattable
        {
            if (rowCount < 1)
            {
                throw new ArgumentOutOfRangeException(nameof(rowCount), Resources.Numerics.ArgumentMustBePositive);
            }

            // Verify Source

            if ((uint)columnIndex >= (uint)ColumnCount)
            {
                throw new ArgumentOutOfRangeException(nameof(columnIndex));
            }

            if ((uint)sourceRowIndex >= (uint)RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(sourceRowIndex));
            }

            if (sourceRowIndex + rowCount > RowCount)
            {
                throw new ArgumentOutOfRangeException(nameof(rowCount));
            }

            // Verify Target

            if ((uint)targetRowIndex >= (uint)target.Length)
            {
                throw new ArgumentOutOfRangeException(nameof(targetRowIndex));
            }

            if (targetRowIndex + rowCount > target.Length)
            {
                throw new ArgumentOutOfRangeException(nameof(rowCount));
            }
        }

        #endregion


        //}}@@@
    }


}
